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1 INTRODUCTION 
1.1 Motivation 



The coalescence of stellar mass black hole-black hole 
(BH-BH) binaries is one of the most anticipated sources 
of gravitational waves for ground-based interferometers 
such a s LIGQj or VIRGCn. Whether formed primor - 
dially l|Belczvnsk i et alj l200gr'2004': 'Sadowski ct al."2008') 
or in dens e star cl uste r s (IPort cgics Zwart & McMillan 2000; 
Giiltekin et al.1 booi: lO'Learv et all l2006l: iGiiltekin et all 



ABSTRACT 

Stellar mass black holes (BHs) are expected to segregate and form a steep density cusp 
around supermassive black holes (SMBHs) in galactic nuclei. We follow the evolution 
of a multi-mass system of BHs and stars by numerically integrating the Fokker-Planck 
energy diffusion equations for a variety of BH mass distributions. We find that the BHs 
"self-segregate" , and that the rarest, most massive BHs dominate the scattering rate 
closest to the SMBH (< 10^^ pc). BH-BH binaries form out of gravitational wave 
emission during BH encounters. We find that the expected rate of BH coalescence 
events detectable by Advanced LIGO is ~ 1 — 10^ yr"-'^, depending on the initial mass 
function of stars in galactic nuclei and the mass of the most massive BHs. We find 
that the actual merger rate is likely ~ 10 times larger than this due to the intrinsic 
scatter of stellar densities in many different galaxies. The BH binaries that form this 
way in galactic nuclei have significant eccentricities as they enter the LIGO band (90% 
with e > 0.9), and are therefore distinguishable from other binaries, which circularise 
before becoming detectable. We also show that eccentric mergers can be detected to 
larger distances and greater BH masses than circular mergers, up to ~ 700 M©. Future 
ground-based gravitational wave observatories will be able to constrain both the mass 
function of BHs and stars in galactic nuclei. 

Key words: galaxies:kinematics and dynamics-galaxies: nuclei-black hole physics- 
gravitational waves 



through the LIGO data stre am looking for such circular in- 
spirals (I Abbott et al.l 200a). but might miss many eccen- 
tric events (IMartel fc Poissonlll999l : lTessmer fc Gopakumarl 
120071 : iMandel et al.ll2008h . In principle, eccentric inspirals 
are well suited for detection as well, however they have not 
been expected to be a signifi cant source for terrestrial de- 
tectors (except see lWenluOOa ). 



20061 : lO'Learv et aLliioOTl : iMiller fc Lauburej|2008l ). nearly 
all of these binaries are expected to be circularised by 
gravitational wave (GW) emission before they are de- 
tectable by gro u nd ba sed observatories IjQ'Leary et al.ll200q : 
iGiiltekin et al.l [200y). In order to be detected with a 
high signal-to-noise ratio, matched filtering algorithms sift 



* E-mail:roleary@cfa. harvard.edu 
f E-mail:bkocsis@ias.edu 
I E-inail:aloeb@cfa. harvard.edu 
^ http://www.ligo.caltech.edu/ 
^ http://www.virgo.infn.it/ 



In this paper, we propose an important additional 
source of gravitational wave sources, which preferentially 
form binaries that are still eccentric as they merge. In 
galactic nuclei with supermassive black holes (SMBH) of 
mass MsMBH < 10^ M0 , relaxation times are often less 
than a Hubble time, and can result in the formation of 
steep density cusp of stars and stellar mass black holes 
(BHs). Indeed, as many as ~ 20,000 BHs are expected to 
have segregated into the in ner ~ 1 pc of the Milky Way 
(Morris 1993: Miralda-Escude fc Gouldll200Gl : iFreitag et ahl 
I2OO6I : iHopman fc Alexandedl2006bl . hereafter HA06). When 
two BHs have a close encounter, they can release sufficient 
amount of energy in GWs to form a tight binary, and merge 
in less than a few hours. These mergers are almost always ec- 
centric, in contrast to circular inspirals expected in globular 
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clusters. A similar process can occur in the runaway growth 
of BHs when no SMBH is present that will also result in 
eccentric mergers; however, such systems are inheren t ly un- 
stable and not long-hved (JQuinlan fc Shapirol Il987l . 1 19891 . 
ll990l : lLe3[l993 ). 

The GWs generated during the evolution of these or- 
bits are endowed with a rich structure. For sufficiently small 
impact parameters, the GWs can be detected with sec- 
ond generation GW instruments already during first passage 
IjKocsis et al.l 120061 ). Subsequently, the eccentricity is close 
to unity for several thousand orbits after the formation of 
the binary, generating a train of short-duration distinct GW 
bursts with a continuous frequency spectrum. As the eccen- 
tricity gradually decreases, the separation between subse- 
quent bursts decreases, and the waveform becomes contin- 
uous in the time-domain and decouples into discrete har- 
monics in frequency. Interestingly, as we show below, the 
eccentricity is non-negligible all the way to coalescence. The 
full waveform, consisting of a correlated set of distinct GW 
bursts, which evolves into a continuous eccentric inspiral sig- 
nal, would not be discovered in the GW data by existing data 
analysis techniques such as GW burst search algorithms, nor 
by small eccentricity inspiral templates, since the contribu- 
tion of individual bursts to the GW power is small and the 
contribution of h igher harmonics is significant throughout 
the evolution (see lYunes et al.ll2008l . for a similar discussion 
for extreme mass ratio GW bursts during encounters with a 
supermassive black hole). 

Eccentric inspirals are more luminous in GWs than 
circular inspirals and extend to higher frequencies, and 
so may be detectable to higher redshifts. Furthermore, 
these events may be capable of detecting intermediate mass 
black hol es using terrestrial instruments (iBr own et al.ll2007l : 
iMandel e t al. 2008). The evolution of the signal during the 
initial phase of eccentric inspirals can be described using a 
small number of parameters, allowing one to construct very 
efficient and sensitive data analysis algorithms to search for 
these signals in the LIGO data. Additionally, these signals 
spend a longer physical time in the LIGO band, possibly re- 
ducing the false alarm rate for a fixed number of templates 
that cover the observational period. The modulation intro- 
duced by the Earth's rotation during such an event can be 
utilised to further improve the measurement accuracy for 
the source position. 

The GW waveforms have been studied for eccen- 
tric inspirals through several different approaches. Post- 
Newtonian (PN) waveforms exist up to 1.5PN order be- 
yond the leading order gravitational radiation effects (oc v^) 
for general mass ratio s with spins (jVasuth fc Maian 120071: 
Maiar fc VasuthlbOOsI '). 2PN without spins JDamour et all 



200i, (oc v\ incl;;3i;g the effects of Eccentricity and ra^ 



diation reaction) and 3PN including eccentricity but with- 
out spins and a lso neglecting radiation reaction (oc u®, 
lArun et al.ll2008l . and references t herein), 5.5PN for ex treme 
mass ratios without spins (oc f^^. lTanaka et al.lll996l ). The 
evolution of a binary is chaotic at (or beyond) the 2PN aj 
proximation if t he component BHs have spins ( Levin i200i 
IWu fc Xiell2007^ and also possibly for nonspinning BHs ap 
proaching the unstable cir cular orbit (jPretorius fc Khuranal 
l2007l : IWashik et al] |2008| . i.e. exhibiting zoom-whirl or- 
bits,). Numerical relativity has mostly focu sed on the cir- 
cular inspiral for comparable mass BHs (jBaker et al.l 120071 : 



iBerti et al.ll2007l : iBovle et al. 2008I. and references therein), 
or extrem e mass ratios ( Martelll2004l : IPretorius fc Khuranal 
2007; Ba bak et al.ll2007r ). Recently, however, the waveforms 
have been evaluat ed for a handful of equal mass eccentric 
initial conditions llSperhake et al.ll2007l : IWashik et ai1l2008l : 



_l|b_pei 



[Hinder et aLll2008a 

In this paper, we perform two separate analyses to two 
different problems. First, in §[2l we determine the multi-mass 
distribution of BHs in galactic nuclei for a variety of mass 
models. In § [31 we determine the formation rate of binaries 
due to gravitational wave capture in a variety of galactic nu- 
clei. In §2] we describe the general features of the generated 
GW waveform, calculate the expected signal-to-noise ratio 
for its detection with second generation terrestrial GW in- 
struments and determine the expected detection rate of such 
mergers. We also show that the eccentricity distribution of 
such sources will distinguish binaries formed in the manner 
outlined in this paper from BH-BH binaries formed dynam- 
ically in massive star clusters. 



2 MASS SEGREGATION 



iBahcall fc Wolj (|l976l . hereafter BW76) were the first to cor- 
rectly analyse the relaxation of a stellar population around 
a central point mass, although their analysis was highly 
idealised. They first derived the Fokker-Planck equations 
for a spherically symmetric distribution of a single-mass 
population of low mass stars around a massive black hole 
(Aft ^ Msmbh). They found that after about one half 
of a relaxation timescale the mass density profile of stars 
reaches a steady-state and forms a power-law cusp with re- 
spect to radius, oc r~°', around the massive central object, 
with a = 7/4. In a second paper, they extended their analy- 
sis to look at a multi-mass system and included effects of the 
loss-cone IjShapiro fc Lightmanlll976r ) that results from the 
disruption of stars by the central black hole (JBahcall fc Wolj 
119771 . hereafter BW77). They found that in a two- mass sys- 
tem the more massive objects segregate from the lower mass 
stars by forming a steeper power-law density profile than the 
stars. This is in stark contrast to the evolution of a two-mass 
cluster of stars without a SMBH, w hich can even tually lead 
to the so-called Spitzer- instability l|Spitzeiill969i ). in which 
the high mass stars decouple from the low mass objects en- 
tirely. 

For a population of stars that are sufficiently old (> 
100 Myr), stellar mass BHs with mass ~ 10 Mq are ex- 
pected to become the most massive objects in the system 
and will begin to segreg ate from the main sequence stars 
due to dynamical friction (|Morrislll993l ) . Looking at our own 
Galactic Cent er and using a ty pical stellar mass function, 
Miralda-Escu de fc Gouldl (|200Gh estimated that ~ 2.5 x 10** 
BHs should have segregated into the inner pc. Indeed, HA06 
followed the work of BW77 by solving the time dependent 
Fokker-Planck equations and confirmed these results. They 
found that ~ 18, 000 BHs of 10 Mq each, should segregate to 
form a steep density cusp around the SMBH in our galaxy, 
Sgr A*, with qbh ~ 2, and that the stars and fighter com- 
pact objects_Jorm a shal lower cusp with a* fa 1.4. 

iFreitag et al.l (|2006l , hereafter FAK06) have so far pro- 
vided the most robust analysis of the segregation and dis- 
tribution of a multi-mass system of BHs and stars in the 
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nucle i of galaxies (for a review, see lArtiaro-Seoane et alj 
l2007l 'l. The authors combined large A^ ~ lO'' Monte- 
Carlo Fo kker-Planck simulations with a population synthe- 
sis code (IBelczvnski et al.1 I2OO2I ) to simulate the segrega- 
tion and subsequent formation of a density cusp for an old 
population of stars. Their simulations were i n remarkably 
good agreement w ith previous analyses (e.g. iMorrisI 119931 : 
iMiralda-Escude fc Gould 2000 . HA06). In their simulations 
of a Milky- Way-like galactic nucleus they found a similar 
overall number of stellar mass BHs (~ 20, 000) within 1 pc 
of a 3.5 X 10'^Mq SMBH. However, they found the BHs 
follow a slightly shallower density profile than HA06, with 
Qbh ~ 1.8. Despite the level of detail in the simulations 
of FAK06, the simulations still have their limitations. Be- 
cause of the large range in scales needed to follow the mass 
segregation, and subsequently the large number of particles 
in that volume, each star cannot be modelled by a single 
particle. The densest, most interesting regions of their sim- 
ulations still suffer from small number statistics, and can't 
reveal the precise expectations for the mass distribution of 
stars and BHs. In order to calculate the rate of gravitational 
wave capture events one needs to apply a high-resolution 
method (like HA06) to multi-mass distribution of BHs for 
calculating the rates and detectability of the events. 

Although the BHs are expected to dominate relaxation 
and dynamical encounters in the inner 0.1 pc around Sgr A*, 
their presence has so far eluded detection. X-ray observa- 
tions of our own galaxy and others may provide the first 
window to this inte r esting stellar population. In our own 
galaxy, iMuno et al.l ((20051) have found an overabundance 
of X-ray sources in the vicinity of Sgr A*, which seems 
to fo llow the underlying distribution of stars (jMuno et al.l 
120061 ). For large distances, the X-ray sources appear to fol- 
low the BH distribution as well (FAK06). In other galax- 
ies, many BHs may interact with the dense accretion disks 
around SMBHs, and in some ci rcumstances be intrinsically 
more luminous than the SMBH (I Deegan fc Navakshinll2007l : 
[Nayakshin & Sunvaov 2007). 

The existence of a nuclear population of BHs may 
also be revealed throu gh dynamical encount ers with ob- 
servable stars. Indeed, [Weinberg et al.l (|2005l ) have shown 
that a decade of observations with the next generation of 
30 m class telescopes may detect the scattering of the BHs 
with the so-called 'S-stars' that are within « 0.04 pc of 
Sgr A*. Microlensing by t he BHs is expected to produce 
a much weaker signatu r e |M iralda-Escude fc Gouldl l2000l : 
[Alexander fc LoeblboOll : Ichaname et al. 20011* ). Strong en- 
counters between the BHs and the same population of S- 
stars may dynamically eject the stars with high enough ve- 
locity to escape the potential of the galaxy (jO'Learv fc Loebl 
120081 ) and po pulate the halo wi t h the observed hyper- 
velocity stars (jBrown et al.l |2005| . [200y). Strong encoun- 
ters are also expected to tidally spin up long-lived stars, 
and so the BH population may be inferred through spec- 
troscopic mecisurements of the spin of lower-mass stars 
IjAlexander fc Kumaij I2OOII ). Alternatively, one may infer 
the presence o f the BHs by looking at density distribution 
of old pulsars (jChaname fc Gouldll20o3 ). 



GW interferometer i/S^ or the Earth-based second gen- 
eration GW instruments are ex pected to detect the inspi- 



ral o f such BHs into SMBH s tGair et al 



2004 



Rubbo et al.l 



20081 ) or IMBHs 



bood: iHopman et al.l l2007l: lYunes et al.l 
(JBrown et al.ll2007l : iMandel et al.ll2008l ). respectively, sev- 
eral times per year. Given the test particle limit of the BH, 
such encounters are expected to be a stri ngent test of Gen- 
eral Relativity in the s trong-field regime (|Collins fc HughesI 
l2004l : lArunet al.ll2006l ). 

An interesting method for detecting a dense population 
of BHs is through the GWs they produce as th ey scatter on 
each other and get captured into a tight binarv. lKocsis et al.l 
(2006. ) estimated the expected detection rates of the first 
hyperbolic passage between such stellar mass black holes in 
globular clusters, and found rates to be typically less than 
0.1 yr^^ for second generation GW instruments. Here, we 
consider the detectability of the much stronger GWs that 
follow once the eccentric binary has formed in galactic nuclei 
in this way and evolves toward coalescence. We explore this 
new method in § [3] and [l] However, in order to completely 
determine this rate we need to know how a multi-mass dis- 
tribution of BHs relaxes around a SMBH in galactic nuclei. 
In the remainder of this section we calculate the steady-state 
distribution of BHs for a variety of models. 



2.1 Fokker-Planck Equations 

In our calculations, we follow the analysis of BW77, who de- 
rived the multi-mass Fokker-Planck equations for a spheri- 
cally symmetric and isotropic distribution of stars around 
a SMBH, /A/(r,u) = fuiE), where E = {l/2)Mv'^ - 
GMMsmbh/t" is the mechanical energy of bound orbits in 
the galactic cusp that are outside the loss cone (see below). 
Much of our analysis follows that of HA06, whose terminol- 
ogy we use throughout our calculations. The main difference 
in our work is that we look at a multi-mass distribution of 
BHs with slightly different initial conditions. 

The SMBH determines the motion and dynamics of 
stars within its radius of influence, Vi = GMsmbh/o"*, 
where Msmbh is the mass of the SMBH, and a^ is the 
one-dimensional velocity dispersion of the underlying stellar 
population in the nucleus near the SMBH. For our calcu- 
lations, we define a dimensionless time variable t = t/tr, 
where tr is the relaxation timescale at the radius of influ- 






3(27r/3,/M,)'^/^ 
327r2G2M2n, InA' 



(1) 



where n* is the number density of stars at r^ for the 
dominate population of stars of mass M*, /?* = M,a1, 
and In A « ln(A^sMBH/Af*) is the Coulomb logarithm. 
We also define a dimensionless unit of energy, x — 
— {Mt/M){E/Pt,), and a dimensionless distribution func- 
tion QAiix) = [(27r/3/M*)^/^nr^]/M(S). In these units, the 
Fokker-Planck equation is (BW77,HA06) 



dgM{x,T) 
dr 



5/2 



^Qm{x)~Rm{x), (2) 



Finally, future GW observations with the space based 
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where 



^ M m T'^ , 



max a;, 2; 



-3/2 






(3) 



is the rate stars of mass Ai flow to energies larger than x, 
and Rm{x) is the rate stars are destroyed by the SMBH. 
We are interested in the distribution of BHs in the inner 
pc and therefore look only at the empty-loss cone regime. 
Thus, the angular momentum averaged destruction rate of 
stars of mass Al with energy x is (HA06) 



Rm{x) 



gM(x) 



Tr(x)ln[Jc(x)/JLc] ' 



(4) 



where Tr is the approximate local relaxation time at radius 

r = ri/{2x), 



m: 



T.m9m{x)AP' 



(5) 



Jc{x) — G'MsMBHO"r^(2a;)~^" is the angular momentum 
of a circular orbit at a radius r — ri/{2x), and Jlc is 
the maximum angular momentum for the object to be 
destroyed or consumed by the SMBH. For stars Jlc = 
GMsMBHO-7^(2(a; + a;td))^ a;^^ corresponds to a star com- 
ing within tidal disruption radius of the SMBH, where Std = 
(MsMBH/M.)~^''^ri/i?. and R* is the radius of the star. BHs 
plunge into the SMBH when J < Jlc = 4GMsmbh/c. 

We numerically integrate Eq. ^ with our initial condi- 
tions given in i) l2.2l until QNiix) reaches a steady state for all 
M. For all of our runs this occurs before r = 1. A compu- 
tationally straightforward simple iterative way to generate 
the steady state solution for Eq. ([2| satisfying -^f- = is to 
follow the time evolution of the distribution from the initial 
conditions. 

After the integration is complete, we calculate the num- 
ber density of stars at radius r. 



nAi(r) =2'K'^/'^n, 



<*(>-) 



AxgM{x){4){r) - x 



,1/2 



(6) 



where cf){r) is the dimensionless specific potential. Until now, 
we have assumed that the stars were on entirely Keplerian 
orbits around the SMBH (i.e. that 0(r) = n/r). This as- 
sumption was necessary in order to simplify the Fokker- 
Planck equations to Eq. ([2} (BW76). However, for our cal- 
culations in § O it is important that the BH density goes to 
zero as r approaches infinity. We therefore calculate cf>{r) and 
nM ir) iteratively after calculating the steady-state distribu- 
tion functions gM{x). We do this by sequentially solving for 
nM{r) (Eq. O and 



r 



M{< r)ri 

MsMBHr 



dr 



(7) 



where M{< r) — X^ai /n^ ^Trr nA/(r)dr is the total mass of 
stars interior to r, thus accounting for the stars' contribu- 
tion to the potential as well as the SMBH's contribution. 
We find that we converge to a single solution for (j>{r) and 
nM{r) before about 6 iterations. This solution to the num- 
ber density and potential is fully consistent with the dis- 
tribution functions gM{x). However, Equations [2] & |3] im- 
plicitly assume that, as r — > 00, the potential approaches 



zero and the number density approaches a constant. This 
assumption is approximately followed as long as the core 
radius of the isothermal sphere re = [9a'^/{'iTvGpo)Y 
(jBinnev fc Tremaind 119871 ) is larger than the radius of in- 
fluence of the SMBH. We also note that the number density 
is the statistical average. In cases of the small numbers of 
BHs, we assume that values we determine will match the 
average over many galaxies. However, we caution that the 
dynamics of such systems will likely behave differently than 
assumed here. 



2.2 Initial Conditions and BH Mass Distributions 

The initial conditions of our models are described entirely 
by the distribution functions of unbound stars gM(x < 0). 
Sufficiently far from the radius of infiuence of the SMBH the 
population of visible stars appear to be similar to an isother- 
mal sphere, with a constant velocity dispersio n and a num- 
ber density distribution t hat scales as oc r~^* (jGenzel et al.l 
I2OO3I : ISchddel et al.|[2007l ). We therefore set g*(a;) = exp(a;) 
for a; < 0, which corresponds to an isothermal sphere of 
stars. We look at only one population of stars of mass 
M* = 1 M0 . This is justified by previous analyses of mass 
segregation that show that the underlying stellar popula- 
tion, being considerably less massive than BHs, tends to 
follow the same distribution function independent of mass 
where the BHs dominate the relaxation (HA06, FAK06). We 
confirm this conclusion in § 12.31 Although there is a signifi- 
cant fraction of more massive stars in nuclei, their lifetimes 
are expected to be less than the relaxation timescale at most 
r and are in such small number that they should not con- 
tribute to the dynamics significantly. 

We assume that the BHs initially follow the distribution 
of the stars, and amount to a total fraction, Cbh = X^a/ ^m, 
of the number density. BH natal kick velocities are ex- 
pected to be much less than t he velocity d ispersions of the 
systems we are inte rested in (|White fc van Paradiig Il996l : 
IWillems et al.l 120051 ). and so we expect the BHs' velocity 
dispersion to be similar to the stars'. Therefore, we set 
QAiix) — Cm exp(a;), yielding a BH number density of CMn^ 
near r^ . This is in contrast to the work of B W77 and H A06a, 
who both assume that the population of stars was already in 
thermal equipartition (i.e. M^m = M*/3*). BW77 addressed 
the distribution of stars in globular clusters, which have half- 
mass relaxation times much shorter than the Hubble time 
and so the entire system should reach equipartition. In con- 
trast, in galactic nuclei the relaxation timescale for large 
radii is usually longer than the age of the universe, and so the 
source population could not reach a complete equilibrium. 
HA06 made a numerical error in how they normalised the 
BH distribution and so had reasonable agreement between 
their work and previous works. However, with the proper 
normalisation, they would have found the number of BHs to 
be (Af/Af,)^" = 10^" times larger than their results indi- 
cated. This is because in the case of a velocity dispersion in- 
dependent of mass, the ratio of number densities of objects is 
described by um /n, = Cm /Ct, whereas for a population in 
thermal equilibrium we get nM/n, = CAfM^''^/(C*M*''^). 
Because the steady-state distribution function is most sen- 
sitive to (7(0) and not the functional form g{x < 0), their 



© 0000 RAS, MNRAS 000, 000-000 



GWs from BHs in Galactic Nuclei 5 



results are similar to ours when one uses the incorrect nor- 
malisation. 

One goal of our analysis is to see if the mass distribution 
of BHs may be determined either through dynamical interac- 
tions with their environment, or through the release of GWs. 
We therefore look at a variety of time-independent BH mass 
functions parametrised as power-laws. The models are deter- 
mined by four basic parameters: the total number fraction 
of BHs, Cbh = X^A/ Cm; the slope of the BH mass function, 
/3, where duM /dM oc M~^; the minimum BH mass, Mmin', 
and, finally, the maximum BH mass in the nucleus, Afmax. 
In our calculations, we approximate the distribution of the 
BHs with A^ = 9 discrete distribution functions with 

"A/,-„i„ + {i+l)AA/ 



a 



-fM-^dM, 



(8) 



Af„i„+iAAf 



where AAf = (Afmax — Mmin)/N, and 7 is normalised such 
that "^^Ci — Cbh- We set the mass for each distribution 
function to be the average mass for that bin. 



M, = 



A-'min + {» + l)AA./ 



il/COM 



-13+1 



dM. 



(9) 



We describe the numbers used above in detail below, and 
outline them in Table [l] 

The number fraction of BHs in galactic nuclei is sensi- 
tive to the initial m ass function (IMF) of hig h mass stars. 
For a Kroupa IMF IjKroupa &: Weidnen 120031 ) we estimate 
Cbh ~ 0.001 by assuming that all stars with mass M > 
20 M0 become BHs. A similar distribution of stars formed 
uniformly throughout time was found to be consistent with 
the K- band luminosity distribution of stars in the galactic 
centre ({Alexander fc Sternbereill999l ). However, more recent 
observations of our Galactic Center suggest the IMF slope 
of high mass stars may be cons iderably shallower than those 
in regular clusters (see, e.g., Navakshin fc Sunvaev _ 20051: 
Navakshin et all l2006l: iPaumard et all I2OO6I: iManess et'all 



20071). iPaumard et al.l (2OO3) found that the stars in an 
apparent disk around Sgr A* are best fit with a flattened 
IMF oc m-'''»f with /3imf ~ 0.85 - 1.35 and a deple- 
tion of low mass stars . This is consistent with the work of 
iNavakshin fc SunvaevI l|2005l ). who constrained the number 
of young low mass objects with X-ray observations, and de- 
termined there to be a deficiency of stars with mass < 3M0. 
Most recently, iManess et al.l (|2007l ) looked at late-type gi- 
ants and concluded that the MF of high mass stars may 
be as shallow as A mf ~ 0.85. A shallow MF is also ex- 
pected theoretically ( Morridl 19931 ) and is consistent with the 
most recent hydrod ynamic simulations of star formation in 
the Galactic Center l| Alexander et al.ll2008h . Considering the 
abundant evidence for an alternative IMF in galactic nuclei, 
we therefore look at a wide range of number fractions, set- 
ting Cbh = 0.001,0.01, and 0.1, which roughly corresponds 
to /3imf = 2.3, 1.5, and 0.8 respectively. 

The mass distribution of BHs is relatively uncon- 
strained. Currently, the best constraints come from the cou- 
ple of tens of X-ray bin aries with dynamically determined 
BH masses (|Remillard fc McClintock 2006). Each measure- 
ment often has considerable uncertainty, but the masses 
seem to span a range of ~ 3 — 18 M© . There is now strong 
dynamical evidence for BHs with even greater masses, up 
to 23-34 Mq, in low metall i city environm ents (Orosz ct al. 
I2OO7I : IPrestwich et al.ll2007l : ISilverman fcFilippenko 2008i) . 



Nevertheless, these observations suffer from a severe ob- 
servational bias; the BHs must be in a close binary to be 
observed. Theoretical es timates for the mass d istribution 
are also highly uncertain. iBelczvnski et al.l (|200J) have used 
their sophisticated population synthesis models to determine 
the expected distribution of BH masses in a variety of envi- 
ronments. Typically, they expect most BHs in high metal- 
licity environments to have a uniform distribution of masses 
between ~ 5 — 10 M0, but do find that significantly more 
massive BHs may form from the merger of the BH with 
a high mass companion star. Unfortunately, these massive 
BHs would not be found in X-ray binaries, unless they were 
introduced into one dynamically. Given these uncertainties, 
we choose to parametrise the mass distribution of BHs as a 
power-law duM/dm oc m~^, where we consider 13 — 2 and 
P — 3. We look at the importance of the highest mass BHs, 
by modelling different upper limits on BH mass, Mmax. We 
also h ave a model based on the work of IBelczvnski et al.l 
(|2004| ). In Model BSR, we use a flat /3 = model of BHs 
with masses between 5 and 10 Mq with a fraction ~ 0.01 of 
BHs with mass 25 M©. This is consistent with their Model 
C2 of solar metallicity stars and based on their Figure 7. 

For the remainder of § [2l we focus our results on Milky 
Way like nuclei, and assume that Msmbh ~ 3.5 x 10^ M© 
(iGhez et al1l2005l : lEisenhauer et al.ll2005h and that the ve- 
locity dispersion of the stars and stellar BHs is a, = 
75kms~^ (HA06). In §[3] we also consider other galaxies that 
have relaxation times short enough to reach steady state in 
a Hubble time with Msmbh ~ 10* — 10^ Mq, in order to 
calculate the overall rate of GWs sources in the universe. 



2.3 Results and Implications 

For all of our models, independent of the BH mass function 
(namely, P, and the mass range of BHs), the most massive 
BHs always become the dominant BH species in the inner 
~ 0.1 pc of the galactic nucleus. In Figure [T] we plot the 
number density of the stars and BHs as a function of radius. 
As has been found previously (BW77, HA06, FAK06), deep 
in the Keplerian potential of the SMBH, the distribution 
functions of the stars and BHs become power-laws of the 
negative specific energy of the objects gm{x) oc x^™- , and 
hence have a power-law density profile oc r"'^ "~''™ . We find 
that throughout all of our simulations, the exponent p™ is 
best fit by a linear relationship between the mass ratio of 
the object, m, and the most massive BH, Mmax, 



Pm =Po 



Mrr 



(10) 



where po ~ 0.5 — 0.6. We have found that po usually has a 
small scatter (~ 20%), depending on Mmax and /3; however, 
given po, Eq. (|10|l is often accurate to < 1%. This relation- 
ship is similar to that found in BW77, who found that po ~ 
0.25 — 0.3, when they looked at tw o different components 
with comparable number densities. [Alexander fc HopmanI 
(20081) attributes the steeper density profiles to "strong" 
mass segregation, where the relaxation of the system is de- 
termined by the many low mass objects. In their calcula- 
tions, Po increases monotonically with Afmax, however the 
rate at which it increases is rather slow after po ~ 0.5. In 
contrast, in our simulations with a mass spectrum of mas- 
sive objects and accounting for the loss-cone, we do not find 
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Table 1. BH Models and Results. Unless otherwise noted, these results apply to a Milky Way like galaxy 
with MsMBH = 3.5 X 10 Mq, ct* = 75kms~ and r^ = GMsmbh/""* = 2.7pc. The columns are, from left 
to right, (1) the model name, (2) the slope of the BH mass function, (3) the minimum BH mass, (4) the 
maximum BH mass, (5) the fraction of stars that are BHs, (6) the merger rate of GW capture binaries per 
model galaxy, and (7) the expected AdLIGO detection rate as calculated in § |4] 



Model 


13 


A/min 


J'^i-nax 


riBu/n* 


Merger Rate per Galaxy 


AdLIGO Detection Rate 






(Mq) 


(M©) 




(yr-i) 


(Caoyr-i) 


A 


2 


5 


10 


0.001 


2.2 X 10^10 


5.8 


A/33 


3 


5 


10 


0.001 


2.0 X 10-10 


5.2 


B 


2 


5 


15 


0.001 


2.8 X 10-10 


16 


B/33 


3 


5 


15 


0.001 


2.5 X 10-10 


13 


B-1 


2 


5 


15 


0.1 


5.3 X 10-9 


220 


B-2 


2 


5 


15 


0.01 


8.3 X 10-10 


40 


Be5" 


2 


5 


15 


0.001 


2.8 X 10-10 


15 


C 


2 


5 


25 


0.001 


3.8 X 10-10 


46 


D 


2 


10 


15 


0.001 


3.4 X 10-10 


23 


E 


2 


10 


25 


0.001 


3.2 X 10-10 


12 


E-1 


2 


10 


25 


0.1 


1.2 X 10"** 


1200 


E-2 


2 


10 


25 


0.01 


1.3 X 10-9 


150 


F-1 


2 


10 


45 


0.1 


1.5 X 10-» 


2700 


bsr£ 





5 


10 


0.001 


2.9 X 10-10 


43 



This is evaluated for a SMBH with Msmbh = 1 x 10'^ Mq and a, = 30kms" 
There is an additional population of 25 Mq BHs with ribh/n, = IQ- ^ 
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Figure 1. The number density of stars and BHs for Model B (left) and E-2 (right). The top dash-dotted line is the number density of 
the stars as a function of radius. The alternating dotted and dash-dotted lines show the number density of the separate mass bins used in 
our calculations with the most massive BHs having the largest number density interior to r ^ 0.3 — 0.6 pc. The solid blue line is the total 
BH number density = '^ j^f n i^j (r) . Near r si Ipc the potential of the stars and BHs is equal to the Keplerian potential of the SMBH. 
The rapid drop in the number of BHs drops rapidly here because of the rigid boundary condition at = 0. In all of our simulations, the 
bin of the most massive BHs dominates the number density of BHs in the inner ~ 0.1 pc. 



that the maximum power-law index po to depend sensitively 
on the number or mass of the highest mass objects. How- 
ever, we have not looked at t he large range of parameters of 
[Alexander fc Hopmanl (|2008h . 

Mass segregation in galactic nuclei ceases when the BHs 



begin to dominate the relaxation process in the inner cusp. 
Typically, this is expected when tiehMIh > n,Ml. One 
therefore may try to detect the presence of a cluster of BHs 
through their interactions with the luminous stars in its 
vicinity. We can estimate the relaxation timescale at a radius 
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r from Eqs. ((TJ & ([5]), in terms of the Keplerian potential of 
the SMBH 



tr[r) = 



3(27rt;g(r))^/^ 
327r2G2 1nAEM^M'^A/(r) 



(11) 



where Vc(;r) = ^/GMsmbhJt is the circular velocity around 
the SMBH at radius r. In Figure[2l we have plotted the relax- 
ation timescale as a function of radius for Models B and E-2. 
The BHs begin to dominate the local relaxation processes at 
r ~ 0.5 pc, consistent with the results of FAK06 and HA06. 
We find that the relaxation timescale in the galactic centre 
is highly sensitive to the assumed mass distribution of the 
BHs, since the lower density of the BHs is easily outweighed 
by their greater mass. By finding a tracer of the relaxation 
time in the inner 0.1 pc of the galactic centre, it may be 
possible to determine if such he avy BHs even exist in our 
galaxy (jChaname fc Gould|[20o3 ). 

As discussed above, a population of massive BHs in 
the galactic centre may be revealed through strong en- 
counters with their neighbouring stars. Because the BHs 
are more massive than the typical star, they tend to 
drive stars out of the nucleus, creating a shallow stel- 
lar density profile (FAK06; HA06; for observations se e 
Gcnzol ct al. 200 3i and more re cently ISchodel et al.ll2007^ . 
Miralda-Escude fc Gouldl l|2000l ) proposed that a population 



of such BHs may be seen through stars on very radial or- 
bits that are ejected by c lose encounters by the BHs. In 
some extreme encounters, lO'Learv fc Loebl l|200a ) showed 
that strong interactions between BHs and stars can lead to 
the stars being ejected from t he Milky Way altogether as the 
observed hypervelocity stars (IBrown et al.ll2005l . 120061 ). The 
results of lO'Learv fc Loea (|2008l ) were sensitive to the as- 
sumed BH mass distribution. Self-segregation of the BHs 
is an important consideration when evaluating the likeli- 
hood of the scenar io, and the velocities of ejected HVSs 
jSesana et aLlbOOTJ ). 

As a more extreme example of a strong encount er, BHs 
and stars may physically collide into each other (j Morris! 
1 19931 . FAK06) producing strong X-ray and UV flares in 
galactic nuclei. Interestingly, we find that the total rate of 
physical collisions between stars and BHs can be compara- 
ble to the rate stars are tidally disrupted by the SMBH, and 
may be visible in variability surveys of galactic nuclei. Such 
encounters may have smaller bolometric luminosities than 
stars disrupted by the SMBH due to their lower Eddington 
luminosity. 



3 BINARY FORMATION, INSPIRAL, AND 
EVENT RATES 

When two compact objects have a close encounter, they can 
emit sufficient energy through GWs that they become bound 
and form a binary. To see the importance of such a process 
in galactic nuclei, we first estimate the rate of binary forma- 
tion for two typical 10 M© BHs. We present more detailed 
calculations later in this section. To first order, two BHs 
with an initial relative velocity w will form a bound binary 
due to the emission of GWs if they come within a distance 
Ki 7.4x 10"^ R0(ui/lOOkms"^)"'*''^ at closest approach (see 
Eq. [18] below). This corresponds to an initial impact param- 
eter (the distance of closest approach if the BHs moved in 



straight lines) of (Eq. [17]) b < 2.4R0(w/lOOkms-i)-''''^ 
We expect that a single BH will undergo such an 
encounter approximately every (yrfc^nBHw)"^ ~ 3.4 x 
10^^(nBH/10®pc~^)-^(u)/100kms"^)"/Vr- The total rate 
of binary formation is then approximately (see, e.g., Eq. |31| 1 
/ n|H'''"&^w47rr^dr, where ■w{r) ~ Vc{r) — y/GMsMsaJr 
is the Keplerian circular velocity at radius r. Integrating 
from r = 10~^ pc to r = 1 pc, yields a rate of 3.4 x 
10~^''yr~MNBH/(2 x 10*)]^ for Nbh BHs in a density cusp 
around the SMBH with wbh oc r~^ . This is in good agree- 
ment with our more detailed calculations presented in § 13.21 

In this section, we calculate the formation rate of such 
gravitational wave capture binaries in greater detail. We first 
calculate the conditions necessary to form gravitational wave 
capture binaries and the properties of their eventual inspi- 
ral in different environments. We then calculate the rate 
of binary formation for the multi-mass, segregated clusters 
of BHs that we modelled in § [2] Finally, we extrapolate 
our results for galactic models, including the distribution 
of SMBHs in the universe and accounting for the statisti- 
cal variance in the number density of stars surrounding the 
SMBHs, in order to estimate the comoving rate density of 
binary formation and inspiral. 

In our calculations, we express our equations with the 
symmetric mass ratio of the encounter r) — {mM)/{m+M)^, 
and total mass of the BHs Mtot- With these variables, the 
reduced mass is ^ = Mtot??. 



3.1 Orbital Evolution 

In the following sections, we shall work with dimension- 
less parameters in units G = c = 1, where the total mass 
Aftot has length and time dimensions. We define the di- 
mensionless pericenter distance as Pp = rp/Mtot, dimen- 
sionless semimajor axis for an eccentric orbit with eccen- 
tricity e as a = a/Mtot = pp/i^- — e). In these units, the 
mean orbital angular frequency from Kepler's law is simply 
c^orb = Aft";t a"^''^ = P '" " 
at pericenter passage is tOp = M^lpp 



time is Aio 



f. 



1 

orb 



(1 — e)^", the angular frequency 
(1-l-e)^''^ the orbital 



3/2 



27i"<i'„.-v, . We define the characteristic 



27ra;p , which 



time of pericenter passage as Atp = f~^ 

satisfies Atp — Atorb for circular encounters, but is much 

smaller for eccentric encounters. 



3.1.1 First Passage, Formation of Binaries 

Two BHs of mass m and M will form a binary if they un- 
dergo a close encounter and release enough energy to be- 
come bound, 5E > riMtotW^ /2, where w — \vm — vm\ 
is the magnitude of relative velocity of the BHs at infin- 
ity, and we ignore the presence of the SMBH during the 
encounter and assume that the interaction is local. Given 
the relativistic nature of such events, and the comparatively 
low velocity disp ersions in nuclei, the e ncounters are always 
nearly parabolic (|Quinlan fc Shapiroll987 : Lee 1993 1. In this 
limit, the amount of energy released du ring the encounter is 



umit, tne amount ot energy released du r 
iPeters fc Mathewslll963l : lTurneiJll977l ) 



5E 



857r 



v'M^i' 



12v/2 



,7/2 



(12) 
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Figure 2. The relaxation timescale as a function of radius for Model B (left) and E-2 (right) plotted in Gyr. The solid line shows the 
relaxation timescale for the entire system. The long dashed and dotted lines show the contribution of the BHs and stars respectively. 
Although the stars contribute to the formation of the cusp and are the most common objects in the system, their low mass precludes 
them from dominating the relaxation in the inner ^ 0.6 pc of the cusp. The break in the relaxation time at ~ 0.0002 pc is due to our 
boundary condition that no stars have x > 10*. 



where r^ is the distance of closest approach, 

Mtot 



fe2 

7,2 2 

b w 
2MZ 



1 ^ Ml, 



1* b'^w'^ 
b^w* 



4M-i, 



(13) 



and b is the impact parameter of the encounter, and in the 
second line we have expanded to first order in w/c. 

The final properties of the system are determined by the 
system's final energy, i?flnai = Mtotrjw'^ /2 + 5E, and angular 
momentum, Lgnai ~ Mtotrjbw + 5L, where 



5L 



G-^MtW 



(14) 



is th e amount of angular momentum lost in GWs (jPetersI 
1 1964 ). For nearly all of the encounters we find \SL\ ^ 
Mtotrjbw and so we set 5L = 0. If the total energy of the 
system is negative, Efinai < 0, then the system will remain 
bound with a semi-major axis. 



ao = 



mLv_ 

2£final 



eccentricity. 



eo 



1 + 2- 



ifo2ii)2 



MLv 



(15) 



(16) 



and pericenter, r-po = ao(l — eo)- 

With the criterion Ennai < for binary formation, the 
maximum impact parameter for two BHs with relative ve- 



locity w to become bound is 



/3407r\ 
bm.. - I 3 I 



1/7 



Mt, 






(17) 



Substituting in Eq. (|13|l . this corresponds to a maximum 
pericenter distance 

2/7 2/7 

r„.max = I 7= 1 Mtot-^77^ (18) 



/857r\ ., 



l-K^)'''(4,)^/^»^«/^ 



We shall demonstrate below that event rates are dominated 
by the central regions of the galactic cusp. For m — M (rj = 
1/4) and w = 1000kms"\ 6max « 2900Aftot, and rpmax = 



47 Mtot- The correction in Eq. (|18|) is clearly negligible for 
nonrelativistic initial conditions. 



3.1.2 Eccentric Inspiral 

After the binary forms with dimensionless pericenter dis- 
tance ppo = rpo/Mtot and ecce ntricity eg, it s orbit will decay 
through the emission of GWs (|Peterslll964l ). 



dpp 
dt 



64 M-lv^ (1 - e 



Pp 



)l/2 
V772 



73 2 ,37 4 , , 

24^ +96^ ) (^') 



de 
IE 



15 p4 (l + e)5/2 



1 ^ 121 2 

1 H e 

304 



^-" (20) 



where we assume pp = rp/Mtot ~ a/[Mtot(l — e)] until the 

last stable orbit. 

We evolve the orbits following IPeterd (|l964l ) for the evo- 
lution of the dimensionless periapsis pp, the time-to-merger. 
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f, as a function of the instantaneous eccentricity, e. Note 
that if starting from a parabohc passage, the initial condi- 
tion for the orbit is determined from exactly one parameter, 
PpO . Dividing Eq. (|19|l with Eq. (|20|l yields a separable differ- 
ential equation for Pp{e). The solution satisfying a parabolic 
encounter initial condition pp — ppo at e = 1 is 



Pp{PpO, 



■ PpO«:p(e) 



where 



'^^ ' * 425 



870 
2299 



ei9 



, ^ 121 2 

1 H e 

304 



870 
2299 



on /I / ^ 



(21) 



(22) 



for which Kp = 1 at e = 1. The orbital evolution equa- 



tions (|21ll & (|22p are valid far from the horizon pp ^ 2. 
Once the last stable orbit (LSO) is reached the evolution is 
no longer quasiperiodic and the binary quickly coalesces. In 
the leading order ratio approximation for infinite mass ratio 
and zero spins, 



Pp(eLSo) = 



6 -I- 2eLso 



clso 



(23) 



which we solve numerically for eLSo(ppo) using Eqs. (I2ip & 
^. 

The time evolution of the eccentricity can be written 
using Eqs. dm & (1201) as 



Img(e) 



3480 
2299 



15 /304y--^^ _i 4 

19 I 425 j ^*°*'' P^' 



/Q 



/ 1 I izi ^z 



(1-.2) 



2^3/2 



-de 



(24) 



where img(e) denotes the time remaining until coalescence 
when the eccentricity is e. Close to coalescence this is ap- 
proximately. 



_5_ /304\ 



3480 
2299 



W(e)~^(^l Mtotry-Vpoe^ ife«l. (25) 

while the tot al merger tim e starting from a highly eccentric 



wnue tne tot al merger tm 
initial orbit l|Peterj|l963 ) 



85 M^^.r, 



(1-e) 



2s,7/2 



if e: 



(26) 



Substituting Eqs. (|15|) & p6|l for the initial separation 
and eccentricity ao and eo into Eq. p6p we get, 



3^3 {bw 



,21/2 



170V85^ M^^/\'i/^ 



(27) 



where we have set -Bflnai = SE, ignoring the initial energy of 
the binary. This is an excellent approximation all the way 
out to 0.996bmax bccause of the strong dependence of tmcrgc 
on the initial velocity w. The merger time reaches its maxi- 
mum at 6 = 6max and simplifies to inGMtot/w^: about two 
times the orbital period of a binary with circular velocity w. 
For the most relevant systems, we find that such bina- 
ries will merge before they have a close encounter with a 
third body. The binary can be disrupted if it encounters a 
third body with a closest approach < 2a before it merges. 
For the maximum impact parameter to form a binary, 6max, 
the separation of the binary that forms is a « GMtot/w'^ ■ 
The typical timescale for an encounter to disrupt the binary 




Figure 3. The eccentricity evolution as a function of time to 
merger. All of the orbits start off near e = 1, and depending 
on the initial periapsis ppo, LSO corresponds to the eccentricity 
glso shown with filled circles as marked. Open circles give blso 
for equidistant intermediate values for ppo- The curves should 
not be extrapolated to e < blso a-s the evolution might be very 
different during the BH merger. 



(12770^ 



= W^I{12-kG'^MI 



where the addi- 



tional factor of 3 comes from gravitational focusing. Be- 
cause of the strong dependence Eq. (|27p has on the im- 
pact parameter (ex 6^^"), for any reasonable disruption 
timescale, a slightly closer impact will always compensate 
and allow the binary to merge before being disrupted. In 
fact, in order for the disruption timescale to be comparable 
to merger timescale, the number density of stars must be 
n > 10^^ pc~^ for w — 100 kms"'^. Because every binary 
results in a rapid merger, the rate of inspiral at any time is 
well approximated by the rate of binary formation. 

Equation (|24|l provides an explicit monotonic relation 
between time and eccentricity. This relation is sensitive to 
the free parameters Mtot, rj, and ppo only through the overall 
scaling MtotV'^ PpO- Therefore eccentricity can be thought 
of as a dimensionless time variable. Figure 131 illustrates the 
conversion between t and e. For the parameters, ppo = 40, 
Mtot = 2OM0, and rj = 0.25, the characteristic timescale is 
several hours for moderate-large eccentricities, and minutes 
for the final small eccentricity inspiral. 



3.2 Event Rates 



b' ) 



Given the cross-section for binary formation (ctcs ~ 7rfe„ 
we can now calculate the expected rate of binary forma- 
tion in a single galactic nucleus, Fign. In order to evaluate 
the total contribution of different BH masses, (m,M), ve- 
locities, {vm,VM), and spatial position, r, within a single 
galactic nucleus, we need to integrate the differential rate 
d^FiGN = crcsWfm(r,Vm)fM{r,VM) over the corresponding 
distributions. Here w = \vm — vm\ is still the magnitude of 
relative velocity between the BHs, and fm and /m are the 
6-dimensional distribution functions of the BHs with mass 
m and M derived in § l2l Generally, this calculation requires 
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us to evaluate an integral over 9 variables, 



IGN 



dr4nr 



dM 



dm 



(28) 



//, 



d Vmd VMfm{r,Vm)fM{r,VM)(TcsW, 
,XM>10,J>Jl,C 



where the integration bounds are set consistent with the 
distribution domains (see !j2]). The multi-dimensional inte- 
gration can be greatly simplified to only three variables us- 
ing a few approximations. First in §2.3, we assumed that 
the the velocity distributions were isotropic. For fixed m 
and M, we switch integration variables from {vm,VM) to 
{vm — VM,Vm + VM), and adopt spherical coordinates. Since 
the integrand only depends on w, we can evaluate the inte- 
grals over the other 5 velocity component^ 



fd\^fS 



VM fm{r, Vm)fM{r, V m) (TcsW : 



nm{r)nM{r) I dw VmM(r, w) CTcsW, 



(29) 



where 4'mM{r,w) is the relative velocity distribution which 
can be expressed as a one-dimensional integral using the 
Dirac-5 function for a power-law dimensionless energy dis- 
tribution profile grnix) ex x^'" (see §[5]for a definition). Note 
that the phase space distribution functions of the BHs and 
stars are well approximated by a power-law Qm oc x^"^ for 
X > 10. 

We have numerically integrated the remaining velocity 
integral over the possible range of w for a variety of differ- 
ent slopes Pm and pm and find that to within < 10% the 
integrand is independent of the shape of the relative velocity 
distribution and only depends on the expected value of the 
relative velocities, and can be expressed with the expected 
value of the individual velocity magnitudes as 



dwijJmMirjW) OcsW ^ 7r6,„axfc(r'), 



(30) 



v^{r) using Eq. (fT7)l . and 
nm,{r) ex r-Pm-3/2 -pj^jg jg ^Q jjg expected, since the encoun- 
ters are practically parabolic such that the initial velocity 
profile is negligible compared to the velocity at periastron. 
Thus, Eq. (123 simplifies to 



3.3 Results and Discussion 

The resulting event rates integrated over a single galactic 
nucleus are presented in the sixth column of Table[T]for all of 
our model nuclei. The estimated rates vary between ~ 10~*- 
10~^'^yr~^ over the various models. We discuss the primary 
sources of uncertainties and other important aspects related 
to the event rates in detail below. 



3. 3. 1 Number of Black Holes 

Overall, we find that the merger rate is most sensitive to 
Cbh = tibh/w. , however, the accretion of BHs by the SMBH 
reduces the naive scaling relation ex C|h. Thus, as the num- 
ber fraction of BHs increases by a factor of 100, the rate of 
mergers only increases by a factor of 20. The merger rate 
depends far less so on Mmax and Afmin, since, for larger 
i^^max the total number of BHs tends to be reduced in the 
inner ~ 0.1 pc of the SMBH. The rate, however, remains rel- 
atively unchanged since the cross-section for binary capture 
increases with the BH mass. 



3.3.2 Pericenter Distance Dependence 

For parabolic encounters, the differential of pericenter dis- 
tances that lead to binary formation, according to Eq. (|13|l 
is to leading order 



drr, 



u/bdb 

Mtot ' 



(33) 



where we have fixed w, the relative velocity at infinity. Since 
dr/d6 oc 6 oc dvp/db, the pericenter distribution that lead 
to binary capture (dF/dr^) is uniform out to a maximum 
pericenter distance rpmax. We can calculate the overall dis- 
tribution of encounters that form binaries by changing the 
order of integration of Eq. (|28|) and leaving it as a function 
of r„. 



dP 



IGN 



drv 



drd«;dmdM 47rr nM{r)nm{r) 

T-WTpmM(r,W), 



(34) 



d'^Fi 



^Vc{r)nm{r)nM{r)r 



dr dm dM 
so that the total rate in one galactic nucleus is 



dr 



dM 



dm 



d-^r 



IGN 



dr dm, dM 



(31) 



(32) 



In practice, we calculate F for all of our simulations by cal- 
culating nm{r) and nM{r) for discreet masses and sum over 
m< M. 



^ Some difficulties arise because of the integration bounds depend 
on the other parameters. However, the integrals can be evaluated 
under the approximation rmin ^ r ^ rmax and ^ i> ^ ^'osc(r), 
where rmin and rmax are the minimum and maximum radii for a 
relaxed population of BHs in the nucleus and Desc (r) is the escape 
velocity at radius r. 



where tpmM (r, w) is the distribution function of relative ve- 
locities given by Eq. (|29|l . The limits of the integration de- 
termine the functional form in three main regimes. For rp < 
rpmax[wniax(rmin)l where u)max(r) = 2veBc{r) = 2y^Vc{r) is 
the maximum relative velocity at radius r, Wcsc is the es- 



cape velocity, and Vp 



:;) is defined in Eq. ([181), Eq. IMll 



is independent of Vp and is integrated over the bounds 

<W < Ulmax(r) 
y'min \ r < rmax- V*^^/ 

r or rp niaxi^^max i^rminj J ^ ^p <^ r^ maxi^Wniax\^max J / tUC lim- 
its of integration are determined by Weq{rp), the inverse of 



Eq. (|T8|I and roq(rp) = 2G'MsMBH/wcq (rp) . In this regime, 
the limits of integration are split into two regions 



< m < TOcq(rp) 

^min ^ r ^ rcql^rpj. 



(36) 
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Figure 4. Normalised probability distribution and integrated probability distribution of initial pericenter distances for binaries in galactic 
nuclei (a) and in Spitzer- unstable massive star clusters without a central massive black hole (b; see ^ l3.3.6l for details). To calculate the 
values in this figure, we have assumed that the merger rate is dominated by the 10 Mq BHs with a near Maxwellian velocity dispersion, 
and in galactic nuclei they have density profiles oc r~^, consistent with the density profile of the most massive objects in the nucleus 
(see i) I2.3| l. The solid line is the differential rate of binary formation dV/drp and the dotted line represents the integrated probability 
r(r > Tp). For galactic nuclei, the probability distribution is perfectly uniform out to pericenter distance rpinax('"'max(-001 pc) ^ 12Mtot, 
after which it drops off quickly. For globular clusters without a central massive BH the distribution is uniform out to rp « 250A'/tot , and 
it begins to drop off as a Gaussian profile for r'p > 540Aftot ■ 



and 



Weq(rp) <W < «;max(r) 



(37) 



Finally, for Vp > rpnmx(2\/2t;c(T'max)) the limits of integra- 
tion are 



< TO < TOcq (r-p) 
^inin ^ ^ ^ ^max • 



(38) 



A normalised probability distribution of pericenter distances 
for encounters in galactic nuclei is plotted in Figure |4^, 
where we approximated ipmMirjW) as a Maxwellian distri- 
bution with variance Vc{r). For comparison the distribution 
function for the Maxwellian core of a star cluster without a 
central massive BH is plotted in Figure |4}3. 

3.3.3 Eccentricity dependence 



We can use Eqs. (|2T|I . (|22|) . and ^^ to follow the secular 
evolution of the binary as it passes through the LIGO band 
and merges when it reaches its last stable orbit. In Figure [S] 
we have plotted the secular evolution of formed binaries in 
10% probability intervals, and plotted the eccentricity af- 
ter approximately every complete orbit with small circles. 
Nearly 90% of all binaries actually form within the LIGO 
band, and « 94% have an eccentricity e > 0.3 as it enters the 
LIGO band (when fp — 10 Hz). This signature is unique to 
an active cluster of BHs in a high velocity dispersion environ- 
ment. For lower velocity dispersions as in star clusters, the 



merger distribution of pericenter distances (dF/drp) stays 
uniform out to ~ rpinax(2 x 50kms~^) « 260Mtot- In con- 
trast to galactic nuclei, we expect only ~ 10% of all binary 
GW sources in star clusters to have eccentricities higher than 
0.3 when they enter the LIGO band, and « 8% to form 
within the LIGO band. Binaries which merge due to 3 or 
4-body inte ractions, are even les s likely to have such a high 
eccentricity (|Gultekin et aLlbOOd l. rO'Leary et akl l|2006l ) cal- 
culated the expected distribution of eccentricities for bina- 
ries that form in dense star clusters (see their Fig. 3). In 
their simulations which had over IQQQ mergers from ran- 
dom encounters, they had no cases that had an eccentricity 
> 0.3 when the binary's GWs entered the LIGO band, and 
only 3 cases (« 10~^) with eccentricity > 0.1. Their cal- 
culations also included secular eff ects which could result in 
mergers with higher eccentricity (|Wenl 120031 ). However, in 
these cases, the binary must have merged after one Kozai 
cycle in a hierarchical triple. Since the Kozai cycle is a dy- 
namical effect, it operates on a timescale much shorter than 
the disruption timescale of the cluster. Thus, any merger 
due to the Kozai effect would indicate an active cluster of 
BHs must still exist. This is in contrast to a delayed merger 
from an ejected binary, which takes of order a Hubble time 
to merge. We discuss other aspects of GW detection in § |4] 
below. 

In most instances, the first binary forming encoun- 
ters occur within the LIGO band. However, the S/N of 
such encounters is too small to be detected for low masses 
(|Kocsis et al.ii200Q ). Therefore, encounters which remain ec- 
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Figure 5. The secular evolution of BH binaries. We have plotted 
in 10% probability intervals the evolution of the binaries as they 
decay (M = m = IOMq; see Fig. |4]for details). The solid line 
is the orbit averaged evolution given by Eqs. I I21I I & II22I I. The 
solid circles denote the completion of approximately one orbital 
period, At^rh- The gap between the solid line and e = 1 is due to 
the finite loss of energy during the initial parabolic encounter. For 
nearly ~ 30% of all binaries, the orbit averaged approximation is 
(visibly) not valid, as can be seen by the large space between 
each orbit. The dotted (red) line denotes where the binaries peak 
harmonic is 10 Hz, the lower limit of the LIGO band. Nearly 90% 
of all binaries that form in galactic nuclei are within this limit 
upon first passage. The dashed (blue) line denotes the eccentricity 
at the last stable orbit (Eq. I23II . 



centric throughout the inspiral, especially near plunge, may 
be the most readily detectable encounters. Therefore we are 
interested in the eccentricity of the binary as it reaches the 
LSO (see Eci.[^. From Eq. (PI)| , and the equation of evolu- 
tion (|Peterslll964l . Eq. I24p we solve for the probability dis- 
tribution of eccentricity at the LSO. We plot the eccentric- 
ity distribution at LSO in Figure |5] for both galactic nuclei 
and globular clusters. For encounters which that to direct 
plunge, our calculation gives an eccentricity greater than 
1. However, for normalisation purposes, we include this in 
our calculations, as they interestingly comprise a significant 
fraction of merger events. 

Until now, nearly all LIGO sources were expected to 
have a negligible eccentricit y as they enter the LIGO band 
(but see lMandel et al.ll2008l . for intermediate mass ratio in- 
spirals in star clusters). The comparatively low eccentricity 
binary formed through few-body encounters or standard bi- 
nary evolution circularise before they enter the LIGO band 
and are detected. Therefore, the detection of eccentric in- 
spirals is a strong test of the formation scenario of nuclear 
binaries, and can conclusively reveal the origin of the BHs. 



Figure 6. The eccentricity distribution of events. Plotted is 
the eccentricity distribution of mergers at the last stable orbit 
(dP/dlneLgo) for one example of a galactic nucleus (solid line) 
and a globular cluster (dotted line). Both lines are normalised so 
that they reach a maximum value of 1. blso ^ 1 corresponds to 
encounters that directly undergo a plunge. See the text and Fig.|4] 
for our assumptions and details of the calculation. 



3.3.4 Radius Dependence inside the Galactic Nucleus 

In Figure [T] we plot the cumulative binary formation rate 
for radii larger than r, r(> r), as well as dF/dlnr. For 
most models, the total differential rate of binary formation 
per logarithmic bin is roughly flat. Thus, each logarithmic 
radius interval contributes equally to the rate. We there- 
fore conclude that the rates determined are rather robust 
to the depletion of BHs very clo se to the SMBH as may 
be caused by resonant relaxation IIRauch fc Tremainelll996l : 
iRauch fc Ingallslll998l : iHopman fc Alexandeij|2006al ) or our 
choice of the innermost radius for BHs. In order for the rate 
to be dominated by mergers at large r, the number density 
of the BHs would have to decrease with an exponent r~" 
where a = p + 3/2 < 3/2 (see § [2]). This is precisely the rea- 
son we accounted for the stars in determining the potential 
in Eq. ©, and did not let the density profile of the BHs and 
stars go to a constant value as in previous analyses. 

We do not expect these tight binaries in their subse- 
quent inspiral phase to have any observable effect dynami- 
cally. Overall, we expect ~ 10 — 10^ such binaries to merge 
over a Hubble time. This presents a much smaller source 
of energy than the SMBH, which accretes ~ 10^ BHs over 
a Hubble time (FAK06). However small the intrinsic rate 
in each galaxy, the cumulative merger rate of many galax- 
ies is large enough to be detected by future ground based 
gravitational wave observatories. 
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3.3.5 Cosmological Merger Rate Density 

Using the Msmbh — c* relationship found bv lTremaine et alj 
lJ2002h for higher mass BHs, 



Msmbh « 1.3 x 10* MQ(cr* /200 km s" 



(39) 



we can extrapolate our results from § 13.1.11 to a range 
of galactic nuclei and determine the overall rate of me rg- 
ers in the Un i verse. Observations bv lBarth et al.l ([2005!) & 
I Greene fc Hoi l|2006l ') have demonstrated that Eq. ((39)) ex- 
tends to active SMBHs with masses as low as 10^ M0 , well 
within the range of our interest. We expect the total mass 
of stars within the radius of influence to be comparable to 
the mass of the SMBH, M(< n) « 2Msmbh. If these stars 
follow a radial density profile of n, oc r~^'^ as expected for 
a relaxed system (see §[2|, then their number density at the 
radius of influence is 



n*(A/s 



1.2 X lO^pc"VlO^M0/Afs 



(40) 



This gives a number density about 45% larger than as- 
sumed in § [2] for the Milky Way and used in our cal- 
culates presented in Table [T] but is well within the ex- 
pected scatter in number densities as we discuss later in 
this section. The actual number density at the radius of in- 
fluence for a BH will depend on both the formation his- 
tory of the galaxy as well as the merger history of the 
SMBH. Because the merger rate is usually greatest at 
the smallest radii, we can can determine the rate in any 
nucleus by scaling Eq. (|32|l and evaluating the integral at 
fmin oc \AMsMBH, the radius where the merger timescale of 
the BH into the SMBH is approximately a Hubble timescale 
(approximately the inner radius of BHs). Coincidentally, this 
radius has the same scaling relation to the outer radius, 
J'max = ri cc x/AfsMBH wherc we use Eq. (|39|l . The total 
merger rate, F is simply proportional to oc n|H(^^U')r'^ eval- 
uated at r = rmin. We can approximate the number density 
of BHs as nBH(r) « CBHn*('"i)(r/ri)~^~^'^. Substituting in 
ri, we get nBH(»'i) oc Af^gjj, independent of the power-law 
distribution of BHs, p. The cross-section of binary capture 
times the relative velocity is -Kh^vj oc u)~^^'^ oc M^gjj , 
given w oc (Msmbh/?")^'^. Combining these dependencies, 
we find the merger rate has a relatively weak dependence 
on the mass of the SMBH, F oc Mg(jgjj. Over two orders 
of magnitude in mass, we expect the rate to change only by 
~ 40 — 60%. To test this scaling relationship we have run 
a simulation with A/smbh = 10^ M© and a* = 30kms~^, 
and found that the rate was in fact comparable to the re- 
lationship found here. Any discrepancy is likely due to the 
slight difference in the capture rate of the BHs and stars 
(seeEq. g]). 

Although the merger rate of BHs is not very sensitive 
to the SMBH mass, it is sensitive to the intrinsic scatter 
of n* for each galaxy, which is used to normalise the entire 
distribution of of the BHs (see Eq. [6]). The expected rate 
of mergers is determined by (nj), the average of the number 
densities squared at the radius of infiuence for the distribu- 
tion of galaxies in the universe. Given (nj) — (n*)^ -f a^, 
where an is the variance of the number density at the radius 
of infiuence of a population of SMBHs, we must rescale our 
results by a factor 



r 



{nl) 
(n.)2 



1 + 



(n,) = 



(41) 



Previous studies on GW event rate estimates neglected cor- 
rections due to cosmic variance, Eq. (|4ip . However it is quite 
plausible that not all galaxies with the same SMBH mass 
have e xactly the same num ber of stars in the central cusp. 
In fact. lMerritt et al.l ((20071) determin ed the relaxation times 
for galaxies in the ACS Virgo Survey ((Cote et al.[[2004( '). and 
found a relatively tight correlation, however, there was still 
significant scatter above the mean by a bout an order of mag - 
nitude. We have taken the results from (Merritt et al.[ ((20071 ) 
(specifically all nuclei with ct, < 140kms~^ in their Fig. 1), 
and determined that a^/(n«)'^ ~ 30. Although there is con- 
siderable uncertainty in the actual value of ^, both observa- 
tionally and theoretically, we expect the merger rate to be 
larger than we have so far calculated by ^ > 10-100, and 
scaled our results by ^30 = C/30. In our calculations, how- 
ever, we use the conservative slope of Eq. (|40|) to determine 
the mean value of n <,, and not the shallower slope found by 
(Merritt et al.[ ((2007( ). which would give rate estimates to be 
an order of magnitude larger than found here. 

We calculate the average cosmological merger rate den- 
sity by convolving the rate per galaxy with the number den- 
sity distribut ion of SMBHs in the univ erse. We extrapolate 
the results of (AUer fc Richstond ((20021) . who found the best 
fit number density distribution of massive SMBHs to be 



dns^ 



dAfsMBH 



Ms^ 



M. 



- Msmbh /Al". 



(42) 



assuming this formula is valid all the way to Msmbh = 
10* M0, where Co = 3.2 x 10~" Mg^ Mpc'^ and M. = 
1.3 X 10* Mq. Finally, we get the cosmological merger rate 
by integration. 



7^ = 



FiGN (Msmbh) -TTT dMsMBH 



dMsMBH 



lO^l Mq 
3^lGNC30Mpc"^ 



(43) 



where Fign is the expected rate of mergers for a single galac- 
tic nucleus of a specific model shown in Table (Tj The nor- 
malisation 3^30 Mpc""^ follows from the distribution given 
in Eq. ((42} , and also accounts for the intrinsic scatter of n* 
for a population of galaxies, Eq. (|41|) . For our fiducial model 
of the MW, Model B, we get a comoving rate density of 
8.4 X 10~^°^3oyr~'^ Mpc"''. Because this integral is nearly 
fiat in the log of A^smbh, it results in a similar rate density 
per logarithmic mass bin and is not very sensitive to the 
limits of the integration. Taking the lower limit of currently 
observed SMBH masses, of 10^ M©, we get a rate ~ 60% 
of what we calc u lated here. However, recent observations 
by (Greene et al.[ ((2008[) have shown that the Msmbh — o" 
relation extends to the smallest SMBHs observed, some of 
which reside in galaxies without a classical bulge. Hence, 
there is still ~ 50% uncertainty in the merger rate due to 
the true function of Fign (Msmbh) discussed above, but the 
total rate density is relatively robust given the other uncer- 
tainties in our calculation, especially ^30 and Cbh. 

For future calculations, we define dTZmhi / Appo as the 
rate density for fixed masses m and M analogous to Eq. (|31|l 
and pi)) . 
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Figure 7. Rate of binary formation as a function of the log of the radius for Models B and E-2. The solid line is the cumulative rate of 
binary formation for radii greater than r. The dotted line is the differential rate distribution per logarithmic bin (dF/dlogr). 



3.3.6 Application to Massive Star Clusters 

BHs in massive star clusters without a central massive 
black hole will also undergo an e poch of m ass segrega- 
tion (|Kulkarni et al.l 1 1993: Sigurdss on fc Herng uist 1993), 
in which the BHs segregate to the cluster core, and effec- 
t ive ly decouple from the stars forming their own subclus- 



ter JPortegies Zwart fc McMillanll20od|: iMerritt et al.ll2004 
lO'Learv et al.ll2006l : iMiller fc Lauburj2008l ). This BH sub- 
cluster will continue to interact only with the BHs until 
a sufficient number of BHs are ejected dynamically, that 
they come back into equilibrium wit h the stars. This oc - 
curs approximately when A^bh < 100 (|Watters et al.ll200(J ). 
For massive clusters, with (uibh) ~ 15kms~^ and nen ~ 
10^ pc~^, the cluster evaporates before about < 10* — lO'' yr. 
Scaling Eq. (|32|l to these parameters, the cluster will have a 
BH-BH merger rate of 



A''bh {riBnWBHO-cs) 

, „ ir,-8 -1 / A'^BH \ 

1.8 X 10 yr X 

^ \wooJ 



__nBH_ 
106 pc- 



WBH 



15kms-^ 



-11/7 



(44) 



during this period of evolution, if the number density of 
BHs is uniform in radius. This is c omparable to, but slig htly 
less than merger the rate found bv lO'Learv et al.l l|200q ') for 
the early evolution of a cluster of BHs due to three-body 
and four-body encounters alone. The detection rate of such 
early mergers depends on the number of young clusters (with 
tags < 10* - 10^ yr) within the detection limit of LIGO. 
Globular clusters, an important source of delayed mergers, 
are too old to still have a BH subcluster. However, the young 
clusters in star burst galaxies would be an excellent source 
if they survive sufficiently lo ng in their hosts to un dergo this 
process of mass segregation (jO'Leary et al.ll2007l ). 



The eccentricity distribution of binary capture merg- 
ers in star clusters is plotted in Figure [S] Overall, the rate 
of mergers in star clusters is dominated by small eccentric- 
ity events, which would be detected as circular inspirals by 
ground-based gravitational wave detectors. However, young 
star clusters may have as many, or even more, eccentric 
mergers than are expected in the nuclei of galaxies if there 
are a sufficient number of mergers in young star clusters. 
In massive star clusters, we expect ~ 10% of all gravita- 
tional wave captures to merge with eccentricities similar to 
those in galactic nuclei. Therefore, the distribution of low ec- 
centricity events will be indicative of the source of BH-BH 
mergers, and may be useful in constraining the distribution 
and evolution of BHs in both galactic nuclei and massive 
star clusters. 



4 DETECTION OF GRAVITATIONAL WAVES 

To determine the expected detection rate of sources, we must 
now calculate the maximum luminosity distance to which 
these inspirals are detectable. In this section, we discuss the 
general properties of the waveform, calculate the maximum 
distance of detection, and add up the total expected detec- 
tion rate for second generation terrestrial GW instruments. 



4.1 General Properties of the Inspiral 

The evolution of the binary and the GW signal can be sep- 
arated into three phases: 

[I] Highly eccentric encounters - train of distinct GW bursts 
in time, broadband signal in frequency. 

[II] Moderate-small eccentricity inspiral - continuous GW 
signal in time, dominated by distinct frequency harmonics. 

[III] Merger and ringdown - short duration peak GW power 
and exponential decay. 
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Figure 8. The relevant timescales that determine the GW wave- 
form: the evolution of the orbital time ( Ato^i, ) and the time du- 
ration of pericenter passage ( Atp ) as functions of time to merger 
(tmerger). The Units are minutes on both axes for (M, rj, ppo) = 
(20, .25, 40) and are different for other values as marked. Similar 
to Fig. [3] the orbits should not be extrapolated beyond the LSO 
for a particular ppo, shown with circles. The inspiral is quasiperi- 
odic if Atorb ^ imerger and Atp <S tmcrgcr- The GW signal is 
burstlike if Atp <C Aiorb a-nd is continuous if Atp ~ Ato^b . 



The distinction between Phases I and II can be under- 
stood by studying the evolution of the relevant timescales 
that determine the GW waveform: the orbital time (Aforb = 
oj-j^), and the time duration of pericenter passage (Aip = 
WjT^)- These are plotted in Figure [8] as a function of time to 
merger fmg. For both Phases I and II, |d Aiorb /dtmg| ^ 1, 
so the orbit is quasiperiodic and evolves gradually due to the 
emission of GWs. The characteristic GW emission timescale 
during each orbit is determined by Atp. Figure |8] shows that 
initially (Phase I), Aip <^ Aiorb, implying that the wave- 
form consists of a train of short Aip duration bursts arriv- 
ing quasiperiodically with separation Aiorb- Later, when the 
burst duration timescale Aip becomes comparable to Aiorb, 
the signal becomes continuous in time domain (Phase II). 
Since the orbital evolution is quasiperiodic, the GW signal 
is approximately a sum of discrete frequency harmonics of 
the orbital frequency, /orb- When the eccentricity becomes 
relatively small e < 0.7, (Phase II) the harmonic decom- 
position is quickly convergent, while during Phase I, it is 
more convenient to work with the continuous limit of the 
frequency spectrum. 

During Phase I and II the equations of motion and the 
GW waveforms can be calculated accurately in various ap- 
proximations (see § 11.11 above) . For Phase II, the expected 
signal to noise ratio o f detection has been calculated by 
iBarack fc Cutlen 1120041). and th e e = 1 parabolic case was 
examined in iKocsis et al.l l|2006l ). In the following, we gen- 
eralise these studies to be applicable to Phases I and II. 

Once the orbital separation reaches the last stable orbit, 
the BHs fall in rapidly and form a common horizon. The last 
stable orbit, which marks the end of Phase II, is determined 
by the initial pericenter distance Ppo for an initially parabolic 
orbit shown by circles in Figure [8] The GW waveform dur- 
ing Phase III involves the calculation of the violently chang- 



ing spacetime and the eventual relaxation into a Kerr BH. 
This requires full numerical simulations of the Einstein equa- 
tions. The detectability of the resulting waveforms have been 
examined f or nonspinning binaries and quasic ircular initial 
conditions (JBaker et al.ll2007l : [Berti et allboOTl ). These stud- 
ies have shown that SjN = 10 can be reached up to a dis- 
tance di = 1 - 6 Gpc for total binary mass Mtot ~ 10- 
2OOM0 for AdLIGO. Eccentric m ergers were considered ver; 
recently by [Hinder et~aLl l|2008al lbh and lWashik et all (1200 
for the case of no initial spins, who showed that the resulting 
GW power is comparable to (or sometimes larger than) the 
power released during quasicircular mergers. Future stud- 
ies should address spin effects during the coalescence, they 
might significantly modify the GW power and waveforms. 

We note that the separation between Phases I-II-III 
is valid only if the initial encounter has a minimum sepa- 
ration, rp, that is much larger than the unstable circular 
orbit, ruco ~ 2-4Mtot, depending on BH spins. Direct cap- 
tures, or orbits outside but repeatedly approaching the un- 
stable circular orbit (the so-called 'zoom- whirl' orbits), are 
qualitatively different. Such encounters have been studied 
in the geodesic approximation a ppropriate for extreme mass 
ratios jGair et al.l I2OO5I I2OO6I : IPretorius fc Khuranal l2007l : 
iLevin fc Perez-Gia |2008|) and using full n umerical simula- 
tion for equal masses ([ Hinder et ahlbOOSaT ). In this case the 
GW spectrum is considerably different and the power is con- 
siderably increased. 

The purpose of this section is to derive the signal-to- 
noise ratio for the quasiperiodic Phases I and II and de- 
termine the maximum range of detection for second gener- 
ation terrestrial GW instruments. We leave the assessment 
of Phase III and the zoom-whirl domain to future studies. 



4.2 Signal to Noise Ratio 

Here we briefly review the general calculation of the signal to 
noise ratio for detecting the GW signal, which we can then 
utilise for the wavefor ms generated by GW captu re events. 
We refer the reader to lFlanagan fc HughesI (| 19981 ) for more 
details. 

In general, the signal to noise ratio of a GW detection 
is defined as 



S^ 

N2 



\h'{f,0)\ 
Snif) 



d/, 



(45) 



where h{f, 6) is the Fourier transform of the GW signal 
weighted by the antenna beam patterns, 6 — {6i} are the 
physical parameters describing the source and the detec- 
tor orientation, Sh{f) is the one-sided noise spectral den- 
sity in units of Hz-^, and /min ^ / ^ /max correspond to 
the frequency band of the instrument, e.g. (/min, /max) ~ 
(10, 10'') Hz for AdLIGCQ. 

The sky position and binary orientation averaged root 
mean square signal-to-noise for a single orthogonal arm in- 
terferometric GW instrument is 






hXf) d/ 
5/S4/) / ' 



(46) 



http://www.ligo.caltech.edu/advLIGO/scripts/summary.shtml 
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where hc{f) is the characteristic isotropic GW amplitude 
defined as 



he = 



2^ = 



TvdL V d/ 




(47) 



where dL{z) is the luminosity distance to a source at a cos- 
mological redshift z, dE/df is the one-sided GW energy 
spectral density on a spherical shell at infinity. The second 
equality corresponds to the stationary phase approximation 
for a quasiperiodic signal sharply peaked at frequency /. In 
this case, E is the total GW power at frequency /, which 
evolves slowly in time according to its time derivative /. 

Equations (|46|) and (|47|) can be generalised for a signal 
co nsisting of discrete harm onics /„, with negligible overlap 
as l|Barack fc Cutledl2004l ) 



iV2 



E 



mill '^Jn'^hyjn) Jn 



and 



ndr 



2En 

T 



(48) 



(49) 



where En is the GW power radiated at frequency /„. For 
quasiperiodic orbits with an intrinsic orbital frequency /orb, 
the observed frequency harmonics at redshift z are given by 



Jn — ^/orb,. 



Jot 



(50) 



4.3 Application to the GW Capture Process 

Let us now turn to the detectability of the GWs starting 
from the initial hyperbolic encounter and ending in the vi- 
olent BH merger. Here we derive a computationally more 
efficient equivalent form of Eq. (|48|l . which can be utilised 

for Phases I a nd II. 

Following iBarack fc Cutlen l|2004h . we calculate the bi- 
nary evolution and the GW signal waveform in the lead- 
ing order approximation of : Peters fc M athews ( 1963) , where 
the interacting masses move on quasi-Newtonian trajecto- 
ries and emit quadrupolar radiation. This approximation is 
adequate in terms of the angular-averaged S/N if the ini- 
tial periapse is well outside the unstable circular orbit, e.g. 
Tp S> r\jco- At smaller initial rp, the Newtonian approxima- 
tion underestimates the GW power and is therefore conser- 
vativelfl To lowest order. 



En 



32 



V M, 



10/3 10/3 



'igir 



(51) 



where g{n, e) determines the relative power of the nth har- 
monic for orb its with eccentricity e, which is given by 
IPeters fc Mathewsi (|1963i 'l in the Newtonian approximation 



° Modulations due to general relativistic pericenter precession 
and relativistic beaming would be very important for the real 
data analysis, but not in terms of the the calculation of event 
rates that depend on the angular-averaged total S/N. 



g{n,e) 



32 



+ (1 



2eJn-i H Jn + 2eJn+i 

n 



e)\Jn- 



2Jn + Jn + 2\ 



4 
-I- — — 7' 

+ 3n2 " 



(52) 



Here Ji{x) is the ith Bessel function, and we have suppressed 
the argument x — ne, i.e. Ji = Jiijie) for each i above. These 
waveforms generally have a maximum at frequency ujp for 
e > 0.5 and have a steep cutoff such that the fractional GW 
power beyond bujp is smaller than 10~^ (|Turnerlll977l l. For 
smaller e, the harmonics beyond the first 2uJp are greatly 
suppressed. The number of harmonics necessary to a preci- 
sion 10""^ can be estimated as 



.(1 + e) 



1/2 



^max — 5 — 5 7- V o /o \5o) 

oJoA (1 - e)'^/2 

e.g. rimax = (10, 40, lO") for e = (0.3, 0.7, 0.99), respectively. 
Note that Eq. (|5ip depends explicitly on e. Therefore 
in order to directly evaluate the S/N integral in Eq. (|48|l 
over /„ one needs to invert the frequency evolution equation 
fn{rpo,eQ,e). Computationally it is much more efficient to 
change the integration variable from /„ to e, 



dfn 



-^dln^ 
2 ao 



18 l + je' + je" de 

IQ 1 _ lMp2 _ 121 p4 p ' 
^^ ^ 304 '^ 304^ ^ 



(54) 



and reverse the order of the sum and the integral 



S^ 



iijMtot^PpO 



95 d| 



<:lSO 



E g{n,e)s{e,eo) de 
n^Shifn) e ' 



n = 2 



(55) 

where glso corresponds to the particular ppo (see Eq. [23]), 
eo is given by Eq. (|16|l . and 



s(e,eo) = 



eo 



p2^3/2 



) 1740 
-^^ (l + eg)(l-e^ 
1 _ 183p2 „ 121 „4 
304 304 



(56) 
Note, that Sh{f) = 00 is assumed outside of fmin ^ f ^ 
/max- The upper limit of the sum in this form can be adjusted 
to the required calculation precision using Eq. (|53p . 

For Phase II, we can rewrite Eq. (|55|l by changing the 
sum over n to a continuous integral over /„. 

/ma. 



S^ 

N2 



■ vMtot,zPpo 



95 



rfi 



CLSO /ir 



g{n, e)s(e, eo) d/ de 
foA,zSh{f) f e ' 



(57) 

where n = ///orb.z, /orb.z = /orb(rp,z, e,eo) given above by 
Eq. (I50j, and rp,^ = Mtot,zPpO = (1 + z)MtotPpO- 

In addition to their numerical advantages, Eqs. (|55|l and 
(|57|) can be used to study the time-frequency evolution of the 
instantaneous S/N accumulation rate as the orbit evolves. 
Figure [9] shows the contribution of the first n ^ 100 har- 
monics to the signal-to-noise ratio for AdLIGO (i.e. before 
evaluating the sum or the integral in Eq. I55|) for total masses 
Mtot = (2O,6O)M0 and initial periapse ppo = (10,20,40). 
The figure illustrates the unique frequency evolution of the 
signal consistent with the expectations described above. Ini- 
tially during Phase I, it is broadband in frequency, and de- 
couples into discrete harmonics at smaller eccentricities dur- 
ing Phase I. The contribution of upper harmonics is nonneg- 
ligible even at LSO, especially if the initial periapse satisfies 
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Figure 9. The frequency evolution of GWs as a function of eccentricity for various total mass and initial periapsis as labelled. The 
shading represents the expected signal-to-noise ratio of the first n = 100 harmonics per logarithmic eccentricity bins for AdLIGO. The 
waveform is described by a broadband spectrum at large e (Repeated Burst Phase) that later separates into distinct harmonics as 
eccentricity decreases until the LSO (Eccentric Inspiral Phase). The signal-to-noise ratio is substantial already at e > 0.7. 



PpO < 40. Note that the maximum frequency of the S/N at 
large eccentricities e > 0.8 in Figure |5] is merely a conse- 
quence of not plotting harmonics beyond n = 100, leading 
to a large underestimate of S/N in Phase II. In this case 
Eq. (|57|l becomes more useful than Eq. (|55p . 

The frequency-independent S/N accumulation rate can 
be obtained as a function of eccentricity (or time using the 
e(t) dependence shown in Fig.[3ll, if evaluating the sum over 
n in Eq. (f55)) or the integral over d/// in Eq. ([57)) . but 
not the integral over eccentricity. Figure [TOl shows the result 
for the same masses and periapse as Fig. (9] including an 
additional extreme case of Mtot = TOOM©. For BH masses 
below Mtot ~ 100, the S/N contribution of the high and 
low eccentricity phases are comparable. The S/N during 
highly eccentric encounters dominate for the typical case 
Mtot ~ 2OM0 for ppo < 20, while small e dominates for 
larger masses Mtot ~ 6OM0 for ppo > 20. It is very inter- 
esting that the eccentric encounters between binaries with 
intermediate masses up to Mtot = TOOM© are detectable to 
Si 1 Gpc for close encounters ppo ^ 10 with AdLIGO, be- 
cause these masses are otherwise totally invisible during a 
circular inspiral. This is explained by the broadband nature 
of the signal during Phase I, leading to nonnegligible power 
leaking into the detectable frequency range of the detector 
/ J5 10 Hz. Such massive BHs were thought to be accessible 
only during the violent merger /ringdown pha se for AdLIGO 
(|Flanagan fc Hughes|[l99i : lBaker et al.ll2007^ . 

The total S/N can be obtained by evaluating both the 
sum and the integrals in Eqs. (|55|) and (|57|) . The result can 
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Figure 10. The square of the signal-to-noise ratio per logarithmic 
eccentricity bin as a function of eccentricity for AdLIGO at 1 Gpc. 
Three sets of curves are shown for total masses with different line 
types as labelled and initial periapse ppo = (10, 20, 40) from right 
to left corresponding to Figure [9] Eccentric encounters between 
massive black holes are detectable up to 7OOM0 with AdLIGO 
that arc totally invisible during a circular inspiral. 



be converted into a maximum distance of detection assuming 
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Figure 11. The total signal-to-noise ratio contours as a function 
of the pericenter separation ppo of the first passage and total 
binary mass, Mtot,z for AdLIGO for an equal-mass binary at 
1 Gpc. Different contours show (S^/Af^)^" (averaged over the 
binary orientation and sky position) increasing in steps of 5 for 
the solid and 1 for the dotted lines. For other mass ratios the SNR 
is reduced by {4:T})^'^ . The limit ppo — > oo corresponds to circular 
inspirals in the AdLIGO band, which result in typically a smaller 
S/N than most likely events with ppo < 40. Note that different 
solid contours correspond to a maximum distance of detection 
increasing in steps of Gpc for a detection threshold S/N = 5. 



a detection threshold, e.g. (S^/N^) = 5^. For Phase II, 
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(58) 
and similarly for Phase I. Note, that the integral depends 
on two parameters Mtot and ppo and is independent of the 
mass ratio rj. 

Figure (fTTj) shows contours of d™'"' for masses 2M0 ^ 
Mtot ^ IOOOM0 and initial periapse 4 ^ ppo ^ 100. For the 
hmit of large ppo, the signal is circularised when it arrives 
in the LIGO band, corresponding to the circular inspiral. 
The figure shows that the eccentric inspirals can be detected 
just to a slightly larger distance than the circular inspirals 
for masses Mtot ^ 2OM0 . The difference becomes more pro- 
nounced at much larger masses. A IOOMq-IOOMq eccen- 
tric inspiral is detectable up to 3 x farther than a standard 
IOM0-IOM0 circular inspiral. For even larger masses, d™^'' 
decreases, but remains nonnegligible up to Mtot ~ 700Mq. 
Currently less is known about the existence of such inter- 
mediate mass black holes in galactic cusps, which makes it 
very difficult to make any theoretical estimates on the ex- 
pected rates of such encounters. It is possible that AdLIGO 
will provide the first direct observational limit on the pop- 
ulation of these objects. We explore this along with other 
estimates of the detection rate in S 14.41 



4.4 Detection Rate Estimates 

We set (^S^/N^) ^ 5 as our detection threshold for a 
single GW instrument, and determine the volume averaged 



maximum luminosity distance for detection for each pair 
of masses, d'^'^^{m,M) from Eq. (|58p . assuming a uniform 



distribution of r„ out to r„ 



Thus, the total detection 



rate for Advanced LIGO is estimated to be 
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max ^'Jmax JV'Jmax 

/ "" dppo / / dmdM / 



d^ 



dKo 
dz 



A/min M 



r7^ 



1 + z dppodmdM 



(59) 



where d'^TZ/dppodmdM is the comoving partial binary for- 
mation rate density between masses (m, M) at initial pe- 
riapse PpO, given by Eqs. pip and (|43p . averaged over 
the distribution of AfsMBH and the number density nor- 
malisation n* (see Eq. [42] and 143 P at redshift z, rigai is 
the comoving number density of galaxies at redshift z, 
Zmax is the maximum detectable distance corresponding to 
d.Lina.x{PpO,'m,M), aud dVco/dz is the comoving volume den- 
sity corresponding to the given cosmology. In practice, how- 
ever, we calculate R discreetly and ignore cosmological ef- 
fects, which are not relevant for the next generation of GW 
instruments. 
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d7?.mA/ 47r , jnaxsS 
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'pO 



(60) 



where dTZmM/dppo is the differential rate density of binary 
formation between mass bins m and M (see Eqs. |3l] and 
I43p . The resulting total detection rates are shown in the last 
column of Table 1. Note that the dependence of di, max on the 
binary parameters (ppo, m, M) lead to an observational bias. 
For example, since d^max is relatively larger for ppo ~ 10, 
the detection rate of these encounters is enhanced relative 
to larger and smaller ppo, even though the intrinsic rate of 
these encounters is independent of Ppo- We have plotted the 
differential detection rate as a function of ppo in Figure [T^] 
for all mergers as well as each mass bin. 

Overall, the most massive BHs in galactic nuclei domi- 
nate the detection rate of mergers for AdLIGO. In Figure fTHl 
we have plotted the distribution of detectable mergers as a 
function of radius for the entire population of BHs, as well 
as for each mass bin. The clear domination of the high mass 
BHs is caused by a combination of three important factors: 
(i) their number density is significantly enhanced by mass 
segregation (§[2]); (ii) the signal of the event is much stronger 
for larger masses (e.g., Eqs. 1551 fc I57p : and (in) the cross- 
section for binary capture is greater for larger masses (e.g., 
Eq.[T3. 

As discussed in § 14.31 we expect the actual inspiral of a 
BH with an IMBH may be revealed by AdLIGO if the event 
is sufficiently eccentric during plunge. Although we cannot 
properly account for this in our analysis, we can attempt to 
compare the rate to what we have done in this work. For a 
large mass ratio, m 3> M, the overall cross-section for form- 
ing binaries increases roughly as h^ ex m^^'^. Compared to 
10 Mq BHs, we expect a 1000 M© BH to have a gravita- 
tional wave capture event ~ 3000 times as often as a single 
BH. However, the number of IMBHs in the region is very 
uncertain. If we take the optimistic number of ~ 10 IMBHs 
in a single galactic nucleus as a steady state distribution 
(jPortegies Zwart et al.ll2006l ). then we expect a comparable 
total number of events to IOMq-IOM© BH-BH inspirals. 



© 0000 RAS, MNRAS 000, 000-000 



GWs from BHs in Galactic Nuclei 19 



Our analysis in §[2]suffers from inaccuracies for the most 
massive and rarest BHs. Equation [2] was derived assuming 
a constant density core (and constant relaxation timescale) 
for large r (BW76), which is clearly violated in most galactic 
nuclei. Typically, the total number of BHs in galactic nuclei 
decreased with Mmax. Despite this decrease in number, we 
likely cannot extrapolate our calculations to higher mass 
BHs such as IMBHs, in order to see what effect they have 
on flattening the density profile of stellar mass BHs. 



5 SUMMARY AND DISCUSSION 

In this paper, we have analysed two separate problems. First 
we determined the multi-mass distribution of BHs in galactic 
nuclei, which we then used to analyse the detection rate and 
merger of gravitational wave capture binaries. We integrated 
the time-dependent Fokker-Planck equations for a variety 
BH mass distributions until they reached a steady state. We 
found, consistent with previous results, that within a relax- 
ation timescale at the radius of influence of the SMBH, the 
BHs and stars form approximately power-law density cusps 
(ex r-^/2-Po'"/""""') within the radius of influence of the 
SMBH. Because the BHs are more massive than the stars, 
they have steeper density profiles (oc r~^) and dominate the 
dynamics and relaxation processes in the inner ~ 0.1 pc near 
the SMBH. 

In such dense population environments, the probabil- 
ity of close flybys between two BHs is nonnegligible. Using 
our results for the steady-state distribution of BHs, we cal- 
culated the expected rate and distribution of GW capture 
binaries in galactic nuclei. We showed that after forming, 
these binaries rapidly inspiral and merge within hours. Un- 
like other sources of merging binaries, the BH binaries in 
galactic nuclei form with a characteristic GW frequency in- 
side the AdLIGO frequency band, and the detectable GW 
signal duration is much longer compared to circular in- 
spirals. In addition, for sufficiently small impact parame- 
ters, the waveforms for such eccentric inspirals are broad- 
band, and can be detected for much larger masses, up to 
Mtot > 700 M0 . Th is exceeds the BH masses previously con - 
sidered detectable (JBrown et al.ll2007l : JMandel et al.ll2008l ). 
and as such it opens a new avenue to probe for the existence 
of a population of intermediate mass black holes in galactic 
nuclei. 

Additionally in § 13.3.61 we also estimated the rate of 
GW captures in young, massive star clusters where BHs de- 
couple from stars to form a subcluster in the centre. We 
found that the rate of GW captures in a single cluster may 
be intrinsically larger than in a galactic nucleus. Typically, 
most binaries circularise before becoming detectable, how- 
ever a significant fraction (~ 10%) may still merge with 
residual eccentricity. Given the vastly different eccentricity 
distribution of mergers in galactic nuclei and massive star 
clusters, the dominate source of eccentric inspirals can be de- 
termined with only a few detections. The total rate of such 
events from massive clusters depends on the relatively un- 
constrained number of young clusters within the detection 
limits of AdLIGO. 

Finally, we analysed the properties of the GW signal 
for the eccentric inspiral of comparable mass BHs, start- 
ing from the initial highly eccentric phase up to the fi- 



nal eccentric inspiral. We found that the maximum lumi- 
nosity distance that such events may be detected using a 
single AdLIGO-type interferometer is 1.2 Gpc for compo- 
nent masses m = A4 ^ 10 M0, and up to 3.3 Gpc for 
m = A4 ^ 50 M0. We then used this to calculate the to- 
tal detection rate of signals, by using a model merger rate 
convolved with a realistic population of SMBHs. We found 
that the most massive BHs dominate the detection rate 
of future ground based gravitational wave detectors with 
~ (10 — 10^)^30 events expected per year, where ^30 is a 
measure of the expectation value for the square of number 
densities near the SMBH (see Eq. l4ip . The overall detection 
rate is sensitive to the number fraction of BHs as well as the 
maximum mass of BHs, and so future observations will be 
able to constrain both the average star formation properties 
and upper mass of BHs in galactic nuclei. 

Overall, the expected AdLIGO detection rate for GW 
capture binaries in galactic nuclei is comparable to estimates 
for other sources of gravitational waves. BH-BH binaries 
that form dynamically in massive star clusters may be de- 
tected by AdLIGO ~ 10 — 10'' times per year, depending 
on the total fraction of star formation that occurs in mas- 
sive, long lived clusters ( Portegies Zwart fc McMillanI 200G|: 



Giiltekin et"ai] 2004: O'Learvet all I2OO6I: iGiiltekin et all 



2006"; O' Learv et al.l 120071 : [Miller fc Lauburd I2OO8Y Esti- 
mates for the population of merging binaries that formed 
in the early evolution of present day globular clusters are 
~ 10 — 100 yr~^ alone. Although globular clusters contain 
only a small fraction of the mass of stars in the universe, 
the merger rate of BHs is sufficiently enhanced by dynam- 
ical interactions that the rate is comparable to those ex- 
pected from standard binary evolution for all ot her stars 
(Bclczvnski ct al."2002'.'2004':'Sadowsk i et al.ir2008l : however 
sec Bclczynski ct al. 2007). Merging neutron stars are also 
a promising source of gravitational waves, whose progen- 
itors have been directly observed. Using constraints from 
the modest populat ion neutron star-neutr on star binaries 
known in our galaxv. lKalogera et al.l (|2004a,bj estimate that 
~ 20 — 600 events per year may be detected by AdLIGO. 
We emphasise that the rates of GW captures calculated in 
this paper are independent of other sources of GWs, and will 
occur in addition to other sources. Given the estimates of 
sources from a variety of environments, gravitational wave 
detectors promise to open a new and interesting window into 
the physics, dynamics, and evolution of compact binaries. 

So far, we have not looked at the interaction be- 
tween the BHs and other compact objects such as neu- 
tron stars (NSs). The NSs are not expected to segregate 
significantly in galactic nuclei, however they are intrinsi- 
cally more common than BHs. To test their importance 
we have performed an additional run of Model E-2 to 
look at the expected detection rate of NS-BH mergers. 
We found that the total detection rate of such events is 
about ~ 1% of the total rate (~ lyr~'), and merits future 
study. Gravitational waves from BH-NS inspirals can pro- 
vide interesting constraints on the equation of state for NSs 
(Faber et al. 2002|), given model sources from num erical sim- 
ulations (iFaber et al.ll2006l : iRantsiou et aLll2008l ). To date, 
such analyses have focused on the circular inspiral of BHs 
with NSs, however recent population synthesis studies sug- 
gest the number of BH-NS binaries that form may be as rar e 
as those expected here, ~ lyr~^ (JBelczvnski et al.l 120071 '). 
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Figure 12. The pericenter distribution for AdLIGO detections. Plotted is the differential detection rate (df?/dlnppo) as a function of 
PpO for Models B (left) and E-2 (right). In both figures, the solid line is differential rate for the all detections, the alternating dotted 
and dashed lines are for each individual mass bin M (J^^^^ d/J^jvf /(2dlnppo)) in order of decreasing mass from top to bottom. Note that 
function plotted here is different than in Fig. |4] 
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Figure 13. The AdLIGO detection rate for Models B and E-2. The solid line is the total integrated AdLIGO detection rate at radii 
larger than r, {R{> r)), the alternating dotted and dashed lines are the contribution to the detection rate of each individual mass bin 
(X]m ^mJV/(> ^0/2) in order of decreasing mass from top to bottom. 



Given that GW capture binaries are predominately eccen- 
tric throughout inspiral, we also expect the newly formed 
BH-NS binaries to be eccentric. However, with BH-NS bi- 
naries, the tidal effects of the BH on the NS can provide an 
additional source of energy dissipation in the binary, which 
can enrich the GW signal, and cause it to deviate from the 
point-mass approximations calculated here. Future numer- 



ical simulations with eccentric encounters must be done in 
order to test the importance of eccentricity on the GW sig- 
nal for BH-NS binaries. 

Another area we neglected is the influence of binarity on 
the merger rate of BHs in galactic nuclei. However, we can 
estimate their overall importance by looking at what frac- 
tion of binaries are hard enough to survive in galactic nuclei. 
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Given the distribution of BH binaries in iBelczvnski et alj 
(|2004[) , only a small fraction, < 0.1% of the BH-BH bina- 
ries are sufficiently hard (with an orbital period < 10 days) 
not to be disrupted by repeated encounters with stars or 
BHs. Interestingly, this fraction is comparable to the frac- 
tion of BHs that merge due GW capture. The properties 
of the merger of such binaries have not been explored yet, 
however, given the evidence that such encounters produce 
few eccentric events in massive star clusters, we expect that 
they will not have many eccentric mergers in galactic nuclei. 
Unfortunately, the methods we have used here are not suit- 
able for looking at higher order A'^-body interactions (such as 
exchanges, or binary-binary scattering); simulations similar 
to FAK06 are more adequate for this type of study. The 
rates we calculated might be enhanced if high-dispersion 
environments of protogalactic cores without a SMBH also 
exist, although such sys t ems are old and not long lived 
JQuinlan fc Shapirolll987l , 1 19891 , Il990l : lLe3ll993l ). 

In order to estimate the detection rate of gravitational 
wave capture binaries in galactic nuclei we had to make as- 
sumptions that introduced uncertainty in our calculations. 
For the sakeof clarity, we wish to summarise these uncertain- 
ties and the order of magnitude effect they may have on the 
actual detection rate for AdLIGO. The largest uncertainty 
in the rate calculation is the distribution of BHs in galac- 
tic nuclei, especially the number fraction of BHs (~ 10^) 
and their mass distribution (~ 10) in the centres of galax- 
ies. Observations of our own galactic centre suggest that the 
star formation may be top-heavy, and place the rate on the 
higher end of that reported here. One effect we neglect in 
our calculations is resonant relaxation, which may reduce 
the rate (~ 10) by depleting the central cusp of BHs. Our 
final rate, however, is in practice only logarithmically sensi- 
tive to the inner edge of the cusp of BHs, and the densest 
cusps that dominate our rate will be least affected by reso- 
nant relaxation. The final rate is also nearly logarithmically 
dependent on the population of SMBHs in the local universe 
(by a factor of at most ~ 10), but is much more sensitive to 
the distribution of densities at the radius of influence of the 
BHs (by a factor of ~ 10-10^). Finally, there is some uncer- 
tainty in the detection rate given our calculation of the S/N 
of the mergers (~ 10) since we do not consider corrections 
for zoom-whirl orbits and conservatively ignore the contri- 
bution of the merger and ringdown waveforms of the binary 
to the signal. 

There are many aspects to the merger of GW capture 
binaries that require future work, which will be aided by 
advanced numerical simulations of the evolution of these bi- 
naries. In our calculations we used only the leading order 
(Newtonian) formula to calculate the evolution of the bi- 
naries and the GW waveforms, and do not account for the 
GW signal during the plunge and merger of the binary, nor 
do we consider the contribution of zoom-whirl orbits. The 
recent breakthrough of numerical simulations finally allows 
the full GW signal to be calculated directly. Eccentric merg- 
ers promise to have a richer signal than that for circular bi- 
naries. Future studies should incorporate the contribution 
of the final coalescence to the GW signal, which would in- 
crease the signal strength generated during the inspiral for 
the same capture event, and has a potential to further in- 
crease the maximum range of detection. 

We also found that a fraction (~ 10 — 20%) of mergers 



will occur very close to, or within the estimated last stable 
orbit of BH binaries. These capture events may be zoom- 
whirl orbits, which dramatica lly increases the GW signa l 
compared to a regular inspiral (IPretorius fc Khuranall2007l ) . 
Such encounters are inherently in the strong GR regime, and 
again require a full numerical treatment of General Relativ- 
ity to be studied. A precise assessment of these sources will 
further increase our estimates of detection rates. Our esti- 
mates for the number of GW capture binaries in galactic 
nuclei (and perhaps in massive star clusters) suggest that 
such events may be common enough to be directly detected, 
and may provide the richest gravitational wave sources in 
the universe for ground based detectors. 
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